# Replication Package

**Paper:** "Robust Inference for Interrupted Time Series Models with Serial Dependent Disturbances"

---

## Software Requirements

- MATLAB R2021b or later
- Required toolboxes: Statistics and Machine Learning Toolbox, Econometrics Toolbox

---

## File Description

### Replication Scripts (run these to reproduce tables)

| File | Description |
|------|-------------|
| `production_table1.m` | Reproduces Table 1 (rejection rates of rescaled t-tests) |
| `production_table2.m` | Reproduces Table 2 (bootstrap vs. rescaled t-test comparison) |
| `empirical_wagner.m`  | Reproduces Table 3 (New Hampshire three-drug cap application) |

### Core Functions (required by all scripts above)

| File | Description |
|------|-------------|
| `run_simulation_func.m` | Core Monte Carlo simulation engine (called by production scripts) |
| `compute_ols_inference.m` | OLS estimation and rescaled t-statistics (Eqs. 18–21) |
| `kernelEstimator.m` | Newey–West HAC long-run variance estimator (Eq. 22) |
| `choose_best_ar_model.m` | AR(p) sieve order selection via AIC/BIC |
| `generate_arma_process_pq.m` | ARMA process generator; supports Gaussian, Student-t, and bootstrap resampling |
| `fix_b_asymptotic_cv.m` | Fixed-b asymptotic critical values (Kiefer & Vogelsang 2005) |

---

## How to Replicate

All scripts should be run from the `Replication/` directory in MATLAB.

### Table 1 — Rejection Rates of Rescaled t-Tests

```matlab
cd Replication
production_table1
```

**Settings:** T ∈ {20, 40, 80, 200, 400, 800}, ρ₁ = θ₁ ∈ {0, 0.3, 0.6}, ξ ∈ {0.75, 0.25}, 100,000 replications, no bootstrap.

**Estimated runtime:** 6–18 hours (depending on hardware).

**Output:** `.mat` files saved to `../Results/production/table1/`.

---

### Table 2 — Bootstrap vs. Rescaled t-Test

```matlab
cd Replication
production_table2
```

**Settings:** T ∈ {20, 40, 80}, ρ₁ = θ₁ ∈ {0.5, 0.8}, ξ ∈ {0.50, 0.75}, 2,000 MC replications × 5,000 bootstrap replications.

**Estimated runtime:** 30–90+ hours (bootstrap-intensive).

**Output:** `.mat` files saved to `../Results/production/table2_run###_YYYYMMDD/`.

---

### Table 3 — Empirical Application (New Hampshire Three-Drug Cap)

```matlab
cd Replication
empirical_wagner
```

**Data requirement:** The empirical data are **not included** in this package, as they belong to the original authors:

> Wagner, A.K., Soumerai, S.B., Zhang, F., and Ross-Degnan, D. (2002), "Segmented regression analysis of interrupted time series studies in medication use research," *Journal of Clinical Pharmacy and Therapeutics* **27**, 299–309.

To run `empirical_wagner.m`, please obtain the data from the original paper or contact the original authors, and place the data file at `../Data/wagner_data_2002.txt`. The script will report a clear error message if the file is not found.

**Settings:** T = 30 (31 months, excluding anticipatory spike at month 20), T₁ = 19, bandwidth S_T = ⌊0.75 T^(1/3)⌋, B = 100,000 bootstrap replications.

**Estimated runtime:** Several hours (100,000 bootstrap replications with T = 30).

---

## Notes on Reproducibility

- All Monte Carlo simulations use fixed seeds for reproducibility. The base seed is `8964`. Production scripts for Table 1 and Table 2 use run-specific seeds derived from the run number.
- Random number generation follows MATLAB's default Mersenne Twister (`rng(8964)`).
- Results may differ slightly across MATLAB versions due to differences in floating-point arithmetic or toolbox implementations.
- The bootstrap in Table 2 used 2,000 MC replications and 5,000 bootstrap replications in the published results. The script uses the same settings.

---

## Data Availability

The MATLAB replication code in this package is freely available. The empirical data (Wagner et al. 2002) are not redistributed; see the note under Table 3 above.

---

